Aquatic Gases

Gases are the inputs and outputs of many important processes in riverine systems. However, their ephemeral nature makes directly measuring them difficult, if not impossible.

The most commonly discussed forms of gases take in aquatic systems is ‘dissolved’. A dissolved gas is entirely intermixed within solution. Note that this doesn’t mean that they are evenly mixed through the entire system. It takes time for dissolved gases to diffuse through water, often creating gradients across the water column.

read.csv('data/SSB_GHG.csv') %>%
  filter(sampledate == '2020-01-09') %>%
  select(water_depth, co2_uM_L, ch4_uM_L) %>%
  rename('Depth (m)' = water_depth,
         'CO2' = co2_uM_L,
         'CH4' = ch4_uM_L) %>%
  pivot_longer(cols = -`Depth (m)`,
               names_to = 'var',
               values_to = 'val') %>%
  ggplot(., aes(y = `Depth (m)`, x = val))+
    geom_point()+
    geom_line()+
    facet_wrap(~var, scales = 'free')+
    labs(x = 'Concentration (uM/L)', title = 'Depth Profiles at South Sparkling Bog')+
    theme_minimal()+
    scale_y_reverse()
Profiles of dissolved methane and caron dioxide concentrations at South Sparkling Bog in Woodruff, WI on 1/9/2020. Note how gradients form with depth due to biological processes.

Profiles of dissolved methane and caron dioxide concentrations at South Sparkling Bog in Woodruff, WI on 1/9/2020. Note how gradients form with depth due to biological processes.

Gorsky, A.L., H.A. Dugan, and N.A. Lottig. 2021. Snow Manipulation Greenhouse Gas Measurements at South Sparkling and Trout Bog 2020-2021 ver 1. Environmental Data Initiative. https://doi.org/10.6073/pasta/661201fc3172c062bdc46010047031d5. Accessed 2022-02-18.

Oxygen

An aquatic gas with a lot of historic data and focused monitoring interest is dissolved oxygen because it is a key indicator of ecosystem health, with low oxygen events potentially causing major disruptions like fish kills. Loading in some data from the South Platte via the USGS, we can see that it cycles with season and time of day.

plot <- readNWISdv(siteNumbers = '06711565',
            parameterCd = c('00010','00300'),
            startDate = '2019-10-01',
            endDate =  '2020-09-30') %>%
  rename('DO (mg/L)' = X_00300_00003,
         'Temp (C)' = X_00010_00003) %>%
  select(Date, `DO (mg/L)`, `Temp (C)` ) %>%
  pivot_longer(cols = -Date, names_to = 'var', values_to = 'val') %>%
  ggplot(., aes(x = Date, y = val))+
    geom_line()+
    facet_wrap(~var, ncol = 1, scales = 'free')+
    labs(y = '', title = '06711565')

ggplotly(plot)
## Warning: `group_by_()` is deprecated as of dplyr 0.7.0.
## Please use `group_by()` instead.
## See vignette('programming') for more help
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_warnings()` to see where this warning was generated.

This is a dynamic graph. Trying zooming in on a single day. Both cycles are driven by two main factors (1) the change in solubility of gases with temperature and (2) seasonal and daily patterns in photosynthesis and respiration.

How exactly streams metabolism carbon is still an emerging field. If this interests you I highly recommend looking at the streamPulse Project and the work of Professor Bob Hall.

Greenhouse gases

Another category of interest is greenhouse gases. The three primary greenhouse gases that inland waters play an important role in are carbon dioxide, methane, and water vapor. Carbon dioxide is the product of respiring organisms and can be dissolved in from the atmosphere. Methane is a much more potent greenhouse gas than carbon dioxide, but exists in much lower concentrations in the atmosphere. In freshwaters, it is commonly the product of methanogenesis, the process by which bacteria reduce carbon (usually in the form of carbon dioxide or acetate) for energy in anoxic environments.

Although lakes and rivers only take up a small portion of the landscape, they act as outsized collectors and reactors of carbon and nutrients.

Flux

Gases are always moving, flowing from high areas of high concentrations to lower concentrations. We describe the net movement of something out of a system as ‘flux’. Remember from our intro unit that you can draw the bounds of a system anywhere. We could calculate flux of anything out (or in) of anything! Flux of calcium out of a watershed (called ‘solute export’), the flux of students out of a classroom, or the flux of gas from the surface of a river.

A (younger) Nick Gubbins taking flux measurements in Trout Creek, WI. These chambers are the same kind linked in the Eurorun example below. Photo credit Luke Loken.

How to estimate gas flux

Flux = k(Concentration_system - Concentration_surroundings)

The lecture videos cover this, but as a reminder: there are three things we need to know to estimate gas flux.

  • The concentration of the gas outside the system
  • The concentration of the gas inside the system
  • The speed at which the gas can move in and out of the system

The first two are pretty straightforward. We can use any number of methods to measure gas concentrations. Depending on the gas, sometimes samples are taken in air-tight bottles, preserved, and returned to a lab to be analyzed. Other gases (like carbon dioxide and methane) preferentially absorb certain wavelengths of light. For those gases, we can use portable analyzers to measure concentrations at high frequency. We will be using data from a Los Gatos Ultra-Portable Greenhouse Gas Analyzer (UGGA) in this homework.

You can watch an over-hyped promotional video for the UGAA here.

Estimating k

Once we know the gradient of the gas between the water and the air, we need to tackle that third item on our list. The speed at which gases can move in or out of water is called the ‘gas transfer velocity’ or k. It is largely a function of how ‘open’, ‘rough’, and ‘turbulent’ the water surface is. I use quotes on those because research shows that k is way more complicated than any (or all) of those things. What we do know is that k varies widely across and within aquatic systems.

The problem in measuring k is a classic one in science: measuring water disturbs water. To overcome this, researchers leverage the fact that net flux can be measured directly, and be used to solve for k.

k = (Concentration_system - Concentration_surroundings)/Flux

Chamber measurements

One way researchers measure flux is with the appropriately named ‘flux chambers’. These were originally developed for soil science, but have recently been applied to water as well. The idea is you constantly sample air from a chamber that is sitting just on top of the water surface, you can observed the change in concentration as an estimate of k. Below is a conceptual diagram of how this works. In this homework, we will be working with a chamber that took replicate measurements.

Example of chamber response for CO2.

This diagram is from Eurorun, a large gas flux project. Before you put the chamber down on the water, the gas inside will have the same concentration as the atmosphere (see the graph starting at about 400 ppm CO2). When you put the chamber down on the water it forms a seal and the gradient begins to equillibrate. The slope of the line in the ‘exponential phase’ is k. As gas fluxes in or out of the water, it will eventually approach equilibrium with the surface water. This is the ‘stationary phase’ and can be used as an estimate for the concentration of the gas in the surface of the water. In this assignment, we won’t be using the stationary phase and instead an measurement will be provided to you.

Sources and sinks

Like I mentioned in my flux introduction, fluxes can be across anything and in any direction. One way we can talk about flux is if something is a ‘source’ or a ‘sink’ of what we are interested in. Sources will have a positive gradient and sinks will have a negative gradient. We only need k to tell us the magnitude of the effect!

Assignment

In this assignment, we are going to look at chamber data from the Yahara River in Madison, WI, see if it is a source or sink of methane in the landscape, and make a flux estimate. The Yahara starts north of Madison, flows through a lot of agricultural land (WI is the Dairy State!), with the gauge itself situated on the east side of downtown Madison. The river is interrupted by beaded, glacial lakes.

The code below will generate a site map, but you will need to do some work to get your document to knit with it included. Ask Matt or I how to enable this.

# Note: need to github install mapview before uncommenting this code!
# tibble(lat = 43.089444, lon = -89.360833, Site = 'Yahara River Gauge') %>%
#   st_as_sf(coords= c('lon', 'lat'), crs = 4326) %>%
#   mapview()

The measurement site is co-located with USGS site 05428500 if you’d like additional context. Here is a recent site timelapse courtesy of the USGS.

Generally when working with gases, it is important to adjust for temperature and pressure. For the sake of this assignment, assume all questions are asking for results in ambient environmental conditions.

Q0:

Just like last time, let’s make sure R knows where we are working. Call ‘here()’ to check.

here()
## [1] "C:/Users/gubbi/Dropbox/My PC (DESKTOP-RMMTJHU)/Documents/WR 418/WR-418-assignments/gas_flux_R"

Q1

Now let’s read in our data. This data file is from a greenhouse gas analyzer and has multiple measurements in it at a single site on the Yahara. The file starts with the chamber hanging in the air.

Part A

Read in the ‘gga_03Jun2014_f0000.txt’ file from the ‘data’ folder. Reduce your data to only chamber methane concentrations and datetime. Do not edit your original data file.

gas <- read.csv('data/gga_03Jun2014_f0000.txt', skip = 1) %>%
  filter(!is.na(X.CH4._ppm)) %>%
  mutate(Time = mdy_hms(Time)) %>%
  select(time = Time, methane_ppm = X.CH4._ppm)

Part B

Plot your data as a time series with appropriate labels.

ggplot(gas, aes(x = time, y = methane_ppm))+
  geom_point()+
  labs(x = 'Time', y = 'CH4 (ppm)')

Part C

How many estimates of k can we make with our data? In other words, how many times did the technician raise and lower the chamber?

Q2

Let’s try to estimate k on our own for the Yahara.

Part A

To estimate k we will fit linear models to the curve right after the chamber has been set on the water surface and the methane is diffusing at a steady rate.

The first step in doing so is excluding data from when the chamber is hanging in the air. This is also useful as a local, atmospheric methane measurement.

Calculate a reasonable ‘baseline’ methane concentration for your data.

(Hint: Baseline values will be low!)

atmo <- quantile(gas$methane_ppm, .15)
atmo
##      15% 
## 1.983725
# or
min(gas$methane_ppm)
## [1] 1.9239

Part B

To estimate k we need to extract the sharp rises from the data, as this is when the chamber is still reaching equilibrium. After a a certain point, the air in the chamber will begin to saturate. This will appear in the time series as fluctuations in methane concentrations. We will also want to filter data while the chamber is in the air.

Filter your data to only the sharp rises, after it has been placed on the water surface, but before the chamber begins to saturate. (Hint: this can either be done through times, values, and/or a new index.)

rises <- gas %>%
  mutate(diff = methane_ppm-lag(methane_ppm)) %>%
  filter(diff > 0,
         methane_ppm > 2.5,
         methane_ppm < 10)

Part C

Plot your new, filtered dataset. Don’t worry about making your plot too pretty, just basic labels is fine.

ggplot(rises, aes(x = time, y = methane_ppm))+
  geom_point()+
  labs(x = 'Time', y = 'CH4 (ppm)')

Part D

Find the slope of each segment using a Sen’s slope test. (Hint: look at the Acid Rain homework for an example of this.)

rises <- rises %>%
  mutate(flag = if_else(time < as.POSIXct('2014-06-03 11:15', tz = 'UTC'), 1,
                        if_else(time > as.POSIXct('2014-06-03 11:25', tz = 'UTC'), 3, 2)))

m1 <- rises %>%
  filter(flag == 1) %>%
  na.omit() 

m2 <- rises %>%
  filter(flag == 2)

m3 <- rises %>%
  filter(flag == 3)

s1 <- sens.slope(m1$methane_ppm)
s2 <- sens.slope(m2$methane_ppm)
s3 <- sens.slope(m3$methane_ppm)
s1
## 
##  Sen's slope
## 
## data:  m1$methane_ppm
## z = 18.404, n = 154, p-value < 2.2e-16
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  0.04578588 0.04642932
## sample estimates:
## Sen's slope 
##   0.0460897
s2
## 
##  Sen's slope
## 
## data:  m2$methane_ppm
## z = 16.806, n = 129, p-value < 2.2e-16
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  0.06095954 0.06371653
## sample estimates:
## Sen's slope 
##  0.06263893
s3
## 
##  Sen's slope
## 
## data:  m3$methane_ppm
## z = 17.957, n = 147, p-value < 2.2e-16
## alternative hypothesis: true z is not equal to 0
## 95 percent confidence interval:
##  0.04810378 0.04925871
## sample estimates:
## Sen's slope 
##   0.0487018

Part E

Report the mean of your k estimates.

mean_k <- mean(s1$estimates[[1]], s2$estimates[[1]], s3$estimates[[1]])
mean_k
## [1] 0.0460897

Q3

Part A

You run a water sample from the surface of the Yahara alongside your chamber measurement and find the methane concentration in the water to be 200 ppm. The chamber is 0.5m tall.Calculate methane flux for your site, report your answer in kg of methane per hour per square meter of water surface.

gradient <- (200-atmo[[1]])

# calc flux in ppm/s
flux = mean_k*(gradient)

# convert to kg/hour*m2
flux*(3600/1)*(1/1000)*(0.5)*(0.001) # sec to hour, mg to kg, height of chamber, liter to cubic 
## [1] 0.01642772

Part B

Does our flux show the Yahara River acting as a source or sink of methane?

Part C

Why is there a difference between the concentration of methane in the air and in the stream? Why isn’t it at equilibrium yet?

Part D

Imagine you are tasked with scaling up your estimate of methane flux for the Yahara site from an hour to an entire water year. Assume you only have the funding to make 3 more trips for measurements. When would you take them? What other public data would you be interested seeing to inform your attempt?

Part E

Imagine you are tasked with scaling up your estimate of methane hourly methane flux of the entire Yahara. Assume you only have the funding to make 3 more trips for measurements. Where would you take them? What other public data would you be interested seeing to inform your attempt?

Bonus

Estimate k for USGS site 06711565 on 5/12/2020. Assume a mean station depth of 0.5m.

Note: the bonus for this assignment is meant to be a more open-ended and a stretch! Don’t worry if you are feeling stuck or less than 100% confident in your methods. Add code chunks and text as needed to show your process and justify your methods. A few hints to help you out. * Look into Odum, 1956 in the ‘pubs’ folder of this assignment * Oxygen saturation != oxygen concentration (look into the function ‘DO.saturation()’) * You will need to look up basic meteorological data * Your value for k will be small and positive

# whle method is as per Odum pg 106
# read data
eng <- readNWISuv(siteNumbers = '06711565',
            parameterCd = c('00010','00300', '00095'),
            startDate = '2020-05-12',
            endDate =  '2020-05-13',
            tz = "America/Denver") %>%
  rename('DO (mg/L)' = X_00300_00000,
         'Temp (C)' = X_00010_00000,
         'Spec Cond (uS/cm)' = X_00095_00000) %>%
  select(dateTime, `DO (mg/L)`, `Temp (C)`, `Spec Cond (uS/cm)`)

# Calc mean salin
salin <- mean(eng$`Spec Cond (uS/cm)`)

# Caclulate DO sat 
eng$`DO (%)` <- DO.saturation(eng$`DO (mg/L)`, eng$`Temp (C)`, elevation.m = 1637, salinity = salin, salinity.units = 'uS')

# plot data to check
ggplot(eng, aes(x = dateTime, y = `DO (%)`))+
  geom_point()+
  geom_hline(yintercept = 1)

# get sunrise/set https://sunrise-sunset.org/us/englewood-co/2020/5
sun_up <- as.POSIXct('2020-05-12 05:17:00')
sun_up_pre <- as.POSIXct('2020-05-12 04:17:00')
sun_down <- as.POSIXct('2020-05-12 20:47:00')
sun_down_post <- as.POSIXct('2020-05-12 21:47:00')

# select data just before rise and after set 
eng_pre_rise <- eng %>%
  filter(dateTime < sun_up,
         dateTime > sun_up_pre)

eng_post_set <- eng %>%
  filter(dateTime > sun_down,
         dateTime < sun_down_post)

# plot data to check
ggplot(eng_pre_rise, aes(x = dateTime, y = `DO (%)`))+
  geom_point()

ggplot(eng_post_set, aes(x = dateTime, y = `DO (%)`))+
  geom_point()

# estimate slope after sun down and before sun up
slope_pre <- sens.slope(eng_pre_rise$`DO (%)`)[[1]]
slope_post <- sens.slope(eng_post_set$`DO (%)`)[[1]]

# find deficits
def_pre <- 1-eng_pre_rise$`DO (%)`[1]
def_post <- 1-eng_post_set$`DO (%)`[1]

# set depth
z = 0.5

# calc k
k_bonus <- z*(slope_pre-slope_post)/(def_pre-def_post)
k_bonus[[1]]
## [1] 0.05891587